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ABSTRACT 

The second order analytical approximation of the mean YORP torque components is given 
as an explicit function of the shape spherical harmonics coefficients for a sufficiently regular 
minor body. The results are based upon a new expression for the insolation function, signif- 
icantly simpler than in previous works. Linearized plane parallel model of the temperature 
distribution derived from the insolation function allows to take into account a nonzero con- 
ductivity. Final expressions for the three average components of the YORP torque related with 
rotation period, obliquity, and precession are given in a form of Legendre series of the cosine 
of obliquity. The series have good numerical properties and can be easily truncated according 
to the degree of Legendre polynomials or associated functions, with first two terms playing 
the principal role. The present version fixes the eiTors discovered in the text that appeared in 
Monthly Notices RAS (388, pp. 297-944). 

Key words: radiation mechanisms: thermal — methods: analytical — celestial mechanics — 
minor planets, asteroids 



1 INTRODUCTION 



The Solar radiation affects rotation of minor bodies in the Sol ar system by two principal mechan isms: direct pressure and thermal re-radiation. 
The former has been shown to be of secondary importance i Nesvo mv & Vokrouhlickv^'2008al). because its average effect is negligible. But 
the significance of the thermal re-radiation, known since the paper of Rubincam ( 2000) un der the name of the YORP (Yarkovsky-O'Keefe- 
Radz ievskii-Paddack) effect, is generally acclaimed and empiri cally demonstrated dLowrv et al.ll2007 ; Taylor et al ] |2007l: iKaasalainen et al 
2007 ). A good account of the related problems can be found in ( Bottke et alj2006h . 



The YORP torques acting on an irregular object are rather difficult to model. Until the last year, only the results of numerical simu- 
lations were available, with the torque values generated by summing the c ontribution s from thousand s of facets of a triangulated surface 
I Vokrouhlickv & Capekll2002l:ICapek & Vokrouhlickvll2004l:IScheeresll2007h . Although ScheeresI J2007h provided a formula for the partially 
averaged YORP torque due to a single small triangular facet in terms of the elliptic integrals, his formulation remains of a semi-analytical 
type. The only known general relations between the body shape and the YORP torque consisted in th e statement about the absence of YORP 
for spherical objec ts, and the 'windmill asymmetry' as the necessary condition for its appearance jRubincanjboOOh . As it was shown by 



Breiter et al 



1 20070 . if the thermal conductivity cannot be neglected, the notion of the 'windmill asymmetry' is wide enough to cover even 
the case of ellipsoids of revolution. 

The first step towards the explanation of the YORP torque values in terms of an arbitrary body shape function was made in the recent 
paper bv lNesvornv & Vokrouhlickv! ( 120070 . The body shape was expressed in terms of spherical harmonics and the general expressions for 
the mean component of the YORP torque responsible for the rotation period variations were derived within the second order approximation, 
i.e involving the products of shape harmonics coefficients. The resulting expressions are cumbersome, involving the products of Wigner 
functions for the dependence on obliquity, and tabulated coefficients of power series in cosine of obliquity does not offer an explanation for 
the notion of the YORP order introduced in the paper. Recently, this line of research has been extended to the YORP torque com ponent 
responsible for variations in obliquity, with a new element - the nonzero conductivity influence ( INesvomv & VokrouhlickvluOOSbl) . but at 



the moment of submission we know only the preliminary version of this article, so our comments about its content have to be rather vague. 

Our present contribution, developed in parallel to the work of Nesvorny and Vokrouhlicky, is based upon similar patterns: shape descrip- 
tion in terms of spherical harmonics, a linear, 'plane parallel' temperature distribution model, and the second order approximation. Yet, we 



* Astronomical Observatory of A. Mickiewicz University, Stoneczna 36, PL 60-286 Poznan, Poland, breiter@amu.edu.pl 

t Astronomical Observatory of A. Mickiewicz University, Sloneczna 36, PL 60-286 Poznan, Poland, han@lab.astro.amu.edu.pl 



2 Breiter and Michalska 



believe that a less complicated formulation of the insolation function, that we succeeded to obtain, leads to more convenient final expressions. 
The Wigner functions play the role of a ladder that becomes unnecessary after one climbs up to the final form of the mean YORP torque 
components, expressible in terms of the Legendre polynomials and associated functions o f order 1. Formulating the mean tor que as Legendre 
series, we introduce the notion of the YORP degree, complementary to the YORP order of lNesvornv & Vokrouhlickvl ( l2007h . directly related 
to the degree of Legendre functions used, and helpful when it comes to truncating the final expressions. 



2 GENERAL FORMULA FOR THE YORP TORQUE 

Consider an infinitesimal, outwards oriented surface element dS of a small body, havin g temperature T. If the Lambert emission model is 
assumed, the thermal radiation of this element gives rise to the force i Bottke et"al]|2006h 



d/ = -^— r^ds, (1) 

3 Vc 

where Et is the surface element emissivity, cr is the Stefan-Boltzmann constant, and Vc stands for the velocity of light. Acting on a surface 
element with the radius vector from the centre of mass r, the force produces an infinitesimal torque dM = rxdf. Assuming the common 
emissivity for all surface elements and integrating over the entire surface of the body, we obtain the net YORP torque 

M = -^^(£TUrxdS). (2) 

3 Vc Js 

The temperature T will be considered a continuous (but not necessarily smooth) function of the longitude and latitude of a surface element, 
and of the Sun position. 

Obviously, if the body is spherical, the vectors r and dS are parallel, so the torque M vanishes identically regardless of the temperature 
distribution function. On the other hand, using elementary identities of the vector calculus, we can rewrite Eq. l|2j, replacing the surface 
integral that involves dS by an integral over the volume of the considered body 

M=l— f [Vx(rV)ldV= ^— r (Vr'*)xrdV. (3) 

3 Vc JyL J 3 I'c Jv 

This form is not more useful for the purpose of the present paper than the usual Eq. but it leads to a second basic fact about the YORP 
torque: M vanishes if the temperature distribution is isotropic, regardless of the body shape. In other words, the mean value of the temperature 
over the surface is irrelevant for the computation of YORP. 



3 BODY SHAPE MODEL AND INSOLATION FUNCTION 
3.1 Surface equation in the body frame 

Let us introduce the reference frame Oxyz with the origin at the centre of mass and the axes aligned with the principal axes of inertia. In this 
'body frame' we can express the object shape in the form of a parametric equation 

/ 

r = a + a ^ ^ 0"'(cos 6) [C;.,„ cos mA + S / ,„ sinra/ij , (4) 
l>l m=0 

being the function of the colatitude 9, longitude A, and of the reference radius a. In further discussion we will always use 

M = cos(9, w = sin 6* = Vl - u^- (5) 

The normalized Legendre functions 0™(m) are defined in Appendix lAll Normalized, dimensionless shape coefiicients C/,„ and 5/,„ will be 
assumed small quantities of the first order. Their small values should allow the treatment by means of the perturbation approach, and they 
should guarantee the convexity of figure, because no shadowing is allowed in our model, as long as the Sun remains above the local tangent 
plane. 

The real form of Eq. Q is instructive, but less advantageous in further transformations. For this reason we pass to the complex form of 
the real- valued equation for r 

I 

r = a + a ^ ^ hm(u. A), (6) 

l>\ m=-l 

expressed in terms of spherical functions K;_,„ (see Appendix lAlt . and of complex shape coefficients 
//,() = Cifl, 

fl,m = (C/.;«-''S/,,„)/2, for 771 >0, (7) 
fl,-,n = (-ir(C,,,„ + 7-5,,,„)/2 = (-ir/*„. 
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Inverting Eqs. ([7} we find 

ClO= flfi, Cij„ = (//,,„ + (-!)'" /;,-„,) = //,,„ + /;*„,' Si^m = i[fl,m-(-^y"fl,-m) = i {fl,m-fl,„), 

in agreement with the properties JAllt of harmonic series. 
For the sake of brevity, we introduce a symbol "P, so that 

r = + 

l>l m=-l 

T and its derivatives will be considered small quantities of the first order. 



(8) 

(9) 

(10) 



3.2 Normal vector approximation 

In this paper we will frequently use a modified set of spherical coordinates with associated unit vectors 

wcosA \ , f -siu/i 

w smA , e^ = — = cosA 

W OA 

u { 



dSr 

du 



-MCOS/i 

-u sin A 
w 



(11) 



forming a right-handed orthonormal basis. Restricting the differentiation to the surface of the unit sphere S, we define a spherical surface 
gradient Vs as 



w oA oil 

Given the shape function ([9}, and setting r = re,-, we can compute the outward normal vector on the body surface 



(12) 



dr dr 
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The normal vector A' given above is relatively easy to handle, but the unit normal vector n derived from il3l has so complicated form, that 
we will use its linear approximation with respect to the shape coefficients 



^ N ^ I d'V ^ „ - „ ,T, 

n w OA oil 



(14) 



The first order approximation is sufficient for our purpose, because it will be multiplied by a first order quantity inside the integral l|2j. 



3.3 Solar position 

Although the shape of an object has been described in the body frame, we will temporarily use another reference frame, where the description 
of the Solar position and of the terminator equation is more convenient. Thus we introduce the Ox'y'z' orbital frame, with the Ox' axis 
pointing towards the Sun, Oz' axis directed along the orbital momentum vector of the minor body, and the Oy' axis completing the orthogonal, 
right-handed triad. 

In terms of the usual 3-1-3 Euler angles and rotation matrices R,, the transformation from the body frame Oxyz to the orbital frame 
Ox'y'z' requires (Fig.[T]l 



r' = R3(i?)R,(e)R3(-n)r. 



(15) 



However, a 3-2-3 rotation sequence is preferred in physics and we will use it in the present discussion. As a matter of fact, the two sets of 
Euler angles are very closely related and 



r' = R3(iA)R2(£)R3(0)'-, 

is fully equivalent to Eq. l ll5t if we se0 



2 ^ 2 



(16) 



(17) 



Assuming the circular, Keplerian motion approximation, we use a constant obliquity angle e, whereas Q reflects the 'diurnal' rotation, and 
the argument of latitude i? is an angle that varies according to the 'yearly' orbital motion. The unit vector towards the Sun Mq and the vector 
normal to orbital plane Kq are given by trivial expression in the orbital frame 



no = (1,0,0)"^ n^ = (0,0,i)^ =e'^. 



(18) 



corrected 
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equinox 



Figure 1. Reference frames and orientation angles used in the paper (corrected). 

Applying the inverse (transposed) rotation matrix ( II6I 1 we transform Wq = e\ to the body frame, obtaining 
Wq = (c cos cost/' - sin sin i/r)e.v + (c sin0 cosi/' + cos0sin(/')ev - cosif/Sr. 
Throughout the text we use the symbols 
c = COSE, s = sine, 



(19) 



(20) 



related with the obliquity angle e. Both c and s are considered constant and the motion of the Sun in the body frame consists of the yearly 
motion described by tf/ and the apparent daily motion reflected in (p. 



3.4 Insolation function 

The so-called insolation function is primarily defined in terms of the unit normal vector of the body surface n and the unit vector directed to 
the Sun as 



£; = (l-A)Omax(O,fJ-W0), 



(21) 



where A designates the Bond albedo and cD = OoOg^ is the ratio of the solar constant Oq ~ 1366 Wm~^ to the square of the orbital semi-major 
axis flo expressed in astronomical units. £ determines the energy flux received by a given surface element of an illuminated body. 
According to Eq. ( 114b . we study an approximate insolation function involving 



n-iiQ^er-iie- —r — — w , 

W'^ OA OA Oil Oil 



(22) 



because is independent on u and A. 

The function max (0, h ■ /jq) is actually equivalent to a product H{h ■ Hq) h itQ, where H(x) is the Heaviside t heta function. The values of 
H(h n p)) are: 1 on the sunny side of the body, on the dark side, and ^ on the terminator. As it was demonstrated by Nesvornv & Vokrouhlickv 
1 20070 , replacing the actual terminator by the terminator of a spherical body leads to an error of the second order. Thus, within the framework 
of the first order approximation, we can replace H(n ■ Aq) by H(er ■ «q). Moreover, thanks to the properties of the Heaviside function, we can 
replace H(y)(dy/dx) by d(H(y)y)/dx, so 



max(0, m-Mq) ~ H{er ■ hQ)er ■ Hq- 
= H-(VsH)-(Vs>I'), 

where 

H = H(er-nQ)er-hQ. 



vv^ OA dA re 



1 as 2 as 

dA dA du du 



(23) 



(24) 



In order to expand S in spherical harmonic series, we benefit from the rotational invariance of the scalar product and consider 3' in the 
orbital frame, where 



S' = Hie'i- ■e'^)e'i. -e'^ = H(w' cos A')w' cos A' 



(25) 
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Using the rule l lAlOt we obtain H' = Y^i^k^'i ii^hkW '^') with the coefficients 

^'ij^ = J '^^ J H(w' cos A')w' cos A' Yl,^(u'a')dA' = J du' J \' cos A' Yij^iu' , A') dA' = 

w'e'j{u')du' = -\ &\(u')&Uu')du' = —, (26) 

1 f^l.l J-l 0-1,1 



expressed in terms of the overlap integrals wj'* described in Appendix |A3l and coefficients H^. defined as 

Hj. = J' cosLe*^dL, (27) 
or, explicitly 

(1 ^ fo'- ^'^1' (28) 
[ \n for A:2 = l, 

with the particular case //q = 2, and //2/(r+l = for all \k\ > 0. 

Note that, apart from / = 1, only even degree / and order k terms will have nonzero coefficients but we will not use this information 
explicitly for the time being. 

The transformation of H' back to the body frame is performed by means of the usual Wigner functions apparatus (see Appendix |A2j. 
So we obtain 

/ / ; 



l>0 k=-lm=-l l>0 m=-l 



Z rr, , ^'"'k 



I>0 m=-l 



Throughout the main body of the article, we will use the convention, that £)^ ^ without an explicit argument always designates D^^^ ^(0, £,(/'), 
and Yi^,„ without an argument is a function of the body frame coordinates Yi ,„(u,A). 

If we substitute H from Eq. i29i into ( |23l l, we can write the insolation function as a sum 

fi a= fi<">+£<'\ (30) 



= -(i-A)oXZZ47T^<.^'.' 



l>Ok=-lm=-l 



(31) 



/ J 

£(1) = -d-A)* J] 2 ^^/^•,tf,.„,(Vsy^;A.)-(v,K,,„,), (32) 
l>0 m=-l j>lk=-j 

where we introduced the explicit form (|TOj of the nonspherical shape function f. 

Using a vector calculus identity for the scalar product of two gradients on the unit sphere 

A,(Fi Fz) = V, • V,(Fi F2) = Fi A.Fj + 2(VsFi ) • (V^Fz) + F2 A^Fi , (33) 

we achieve the transformation 

(v.sKj;a) -(V.y/.^) = -i [Yj,kA.Yi,„ + Y/^^A^Yj^k- A,iYj,kyi.m)] = \ [(/(/+ l) + j(7+ 1))K),* Y,,„ + A,{Y j^kYl,m)\, (34) 

where we used the fact, that a spherical harmonic Yi> „y is the eigenfunction of the Laplacian operator on the unit sphere As, with the 
eigenvalue -/'(/' + 1). Now, if we convert the products Yji;Yi ,„ into Clebsch-Gordan series ( IA33K we obtain 

(vsi';,^)-(vsy/,,„) = i^(-i)m(i+ i)+j(j+i)-p(p+ i))g';''^pYpq, (35) 

where the Gaunt coefficients ^'"j*^ are nonzero only if q = m + k, according to their properties {ASOjl. 
Thus the first order approximation of the insolation function is given by 

£«ZZ(fiS-4'i)>''- (36) 

l>0 m=-l 

where 

fiS = -(1-^)*Z^<. = -(1-^)'I'Z^<.^^^^^^^^^ (37) 
= (-l)'"(l-A)02Z/.*Z4^;:;f""Z^<-..^*'^^^^ (38) 



6 Breiter and Michalska 



SpJ = 5 (P(P+1) + J(,/+ !)-/(/+ D). (39) 

Unless otherwise stated, a Wigner d-function without an exphcit argument is d}^^^{e) defined in Eq. JAlSt . Both Eqs. l l37t and ( I38t 
contain many terms that actually vanish in virtue of the properties possessed by functions Hk, , and Q"^ j^^'^' '", but we delay the task of 
purging the series until the last stage of the YORP torque derivation, preferring rather a compact notation in intermediate transformations. 



4 TEMPERATURE MODEL 



If the body is treated as a perfectly nonconducting matter, with the thermal conductivity A" = 0, the temperature of a flat surface element is 
described by 



£tcr 



(40) 



where 5 is the insolation function. This ideal ized temperature function leads to the models of YORP called Rubincam's approximation 
I RubincamlboodlVokrouhlickv & CapeSbooA . 

A commo n practice in modeling the YORP torques for bodies with nonzero coductivity, is to assume a so-called 'pla ne-parallel' tem- 



perature model dCapek & Vokrouhlickvll2004IScheeres & MirrahimilboOSljNesvomv & Vokrouhli"ckil2008bl : lMvsenll2008l) . It is based on a 
one-dimensional (radial) heat diffusion equation 



dT _ d^T 
It ~^'d^' 

where ( measures the depth (( = on the surface, ( > inside the body). The coefficient x is the ratio of the thermal conductivity K to the 
product of the density p and the specific heat capacity Cp 

K 

The boundary conditions on the body surface are related with the conservation of energy, which implies 
dT 



X 



(41) 



(42) 



-fi = 0. 



(43) 



(=0 



The first term is the re-radiated energy, proportional to the fourth power of T, and the product of the Boltzmann constant cr and the surface 
emissivity Sf The second term is responsible for the heat transport towards the body's centre, and the last term is the insolation function fi, 
providing the incoming radiation energy. Moreover, we search only the steady state solution that is a quasi-periodic function of time, and 
such that lim^_>oo ^ = 0- 



Nonlinear boundary conditions are linearized around the constant value Tq, given by 



^4 r-(O) 

e^o-Tq =6] 



(l-A)O 



(44) 



00 4 
leading to the linear model 
T xTo + f, 
where T is given by 
, - \df 

A more detailed description of this model can be found in Bertotti et al. 20031) or I Vokrouhlickv 1998 ). For our purpose, we simply 
state an operational rule, that each term of the insolation function £ having a form 

S„„(M,i)e'(''^+«^\ 



= 8-8] 



(0) 
'00- 



(45) 



(46) 



creates a term fpq defined by 



AeiCT tI fpq = Bpqiu, A) Kpq e' 



Hpip+qip+Vn) 



(47) 



where, assuming a constant rotation rate a> = -0, and a circular orbit with mean motion nj = t//, the thermal conductivity effects are described 
in terms of 

1 



pg 



^l+2Fpq+2Fjq 

sinifipg = sgnipoj -qn^)FpqKpg, 



(48) 
(49) 
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cos 



pp., = -7^^^l " (51) 



= {l+Fpg)Kpg, (50) 

Note, that K^p^^q = Kp^q, but (p-p-q = -fpg- If the second subscript is 0, we will use an abbreviated form 

Kpo = Kp, <ppo = (fip, (52) 

with important special cases K(i= \, and i/jq = 0. 

If we apply the rule l |47| ( to the insolation function series ( I36I38I ) and split t into the spherical part T^^^ and non-spherical correction 
r^^-*, we find 

s,aT'.s,aTt + As,aTl[T^KT^'^)=^^^+Y. Z C + Z Z ^?l^U'n(.u,X), (53) 

l>\ m=-l ' l>0 m=-l 

with 

/ 1.* u 

'fi = -(i-A)o2^i^-— /f„a<,^^"'^'"^^''^'-^'"''\ (54) 
'+> j w\'''Ha 

= ^'-^^'^I.I.I.I.hd-'rs'p.e"-]^^^^^^^ (55) 

pip=Og=-pk=-j 

This temperature model is accurate up to the first power of the shape coefficients fjj^ not only because of the approximate insolation function 
used to derive it, but also because of Tq defined according to the spherical model. 



5 TORQUE INTEGRATION OVER THE SURFACE 
5.1 Approximate integrand 

The integral ^ involves the factor (r x dS). In terms of the modified spherical coordinates in the body frame, 

dS = NdudA, (56) 
where the normal vector A' is given by Eq. ([T3}. Substituting Eq. ([9}, we find 

rxN = a\\+'¥ferX[(l+ f) - Vs'P] = -a^ ( 1 + T)^ x VsT = ( 1 + '¥fL''¥, (57) 

where L'^ = -iiy x Vj is the angular momentum operator i JacksonI 1967 ; Biedenharn & Louck 1981 ) reduced to the unit sphere. Up to the 
second order, we can approximate r x A' as 

rxN~ -ia^ (1 + 2T)L'4' = +21')F. 



According to the action of V' on spherical harmonics in the complex, body frame-fixed basis 



e- = (e.v - ie^) , e+ = e*_ = ^{ejc + /e,,) , 
the vector F = -iL'^'i', is a sum 



1>1 m=-l 



where the vector coefficients 

Fl,m = -' {(rlmhm+ie- + o-linyi,m-ie+ + inYii„eo) , 
involve auxiliary symbols 



{l-m)(l + m+ 1) 
O"/,," = ^1 = 



(l + m)(l-m+ 1) 



1,-m 



(58) 



(59) 



(60) 



(61) 



(62) 



2 l,-itr 1,111 Y 2 

Note that from the point of view of a usual scalar product e±-eo = 0, and So-Sq = I, but e+-e+ = e_ • e_ = 0, and e_ ■ e+ = 1 . 

According to the conclusions drawn from Eq. we can drop Z = m = terms from our temperature model J53b . because they have no 
effect on the YORP torque. Hence, we will use 



(63) 



l>\ m=-l 
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in further considerations. Thus we have reached the stage, where the problem of the YORP torque is reduced to an integral over the surface 
of the unit sphere 

^-1? I^'C^ i [C^-2.|>F.4>)F]n,„d., (64) 
ivc J -I Jo 

where the first term in the square bracket is a first order YORP torque and remaining two define the second order torque. Any second order 
(2) 

term fj^^ in temperature function, depending on the products of shape coefficients, would result in a third order contribution after being 
multiplied by the first order vector F. This justifies the use of the first order approximation of in the second order YORP torque solution. 



5.2 Integrated torque in the body frame 

The first order part of Eq. l |64t is 

^'' = -'^Z t Ztf'^'"'^l I'd. r 7,,n„,d.. (65) 

" l>\ m=-l p\k=-j -^-1 

If we substitute Yjj^ = (-1)*K*^ and use the orthogonality of spherical functions JA8I I. we easily find 
2ia — — 



= _£!^ y y r^- m + tf^ , + mtf^ So] . (66) 

^ V / i / 1 ^ l,m L l,m l,m+ 1 ^ l,m /,m- 1 /,m " J ^ ^ 

;>!;«=-/ 

In the second order we have M^' = M^"^' + with 



(67) 



l>l m=-l p \ k=-j 

resembling Af*'', and 



^(02) = 2 ^ z z z z fp--^^ \\ C ^".^ p'''- (68) 

" m m=-l pi k=-jp>\q=-p 

Of course, the final form of M*^"^ is similar to Eq. ( 166b . with all f^j, replaced by 
^(11) y y \a7j'\,e,+at„P ,^-+mr<"eol, (69) 

whereas M*^"^' looks different, but still the integral of three spherical functions is directly given in terms of the Gaunt integrals ( IA28b . In our 
case, we obtain g'I'f, and &';'f-\ that are nonzero only for m + q + k equal to 0, -1 and 1 respectively (see Appendix I A3 1 . This 

property helps to remove one sum, and 

t ■ 3 I P 

' = > > > fimfDa\cr,„,G, r. ,e-+aJ,„G, ^ t.' ,e++mG, '. ' t. ei)\. (70) 

3v ^—l ^—l ^—l •"•'"•'''•q\ l^pj j,-m-q-l l,m^ l^pj j-m-q+l ^ hpj J,-"i-q j 

l,j.p>lm=-lq=-p 

Further restrictions of summation range, implied by properties l |A30t , will be introduced later. 

The equations derived in this section provide a tool to compute the complete second-order approximation of the YORP torque 

M = M<'^+y[(M*"2^+M*"^ (71) 

if only the spherical harmonics expansion of the temperature distribution function is known up to the first order in shape coefficients. We 
introduce a symbol = 1 as a marker that will help us to locate the contribution of the Af*^"^', neglected by Nesvomv & Vokrouhlickv, c2008bl) . 
in the final expressions. 



6 ELEMENTS OF DYNAMICS 

Given an object on a circular orbit, rotating around its principal axis of the maximum inertia, we can describe its dynamics under the action of 
a torque M in terms of four differential equations related to the 3-1-3 Euler angle^ Q. = -(I>-7t/2, s,& = i//-n/2, and to the angular rotation 



^ corrected 
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rate 



M-ex 



(72) 
(73) 



Meo 

§ = It, = i, (74) 

cM-eo 

a = -<p = oj = ai + c(i?-;js), (75) 

where ^3 is the maximum moment of inertia, is the orbital mean motion and, as usually, c = cose, ^ = sine. The three unit vectors in 
the above equations define an inertial reference frame attached to the body's equator-equinox system]^ ei ~ opposite to the projection of ei 
(normal to the orbital plane) on the equatorial plane Oxy, e2 - opposite to the vernal equinox (i.e. directed along ei x e^), and 63 - directed to 
the north pole. The components of this right-handed, orthonormal basis in the body frame are given by 

ei = -cosfpS x - sim^e,, = - ^ |^e""^e+ + e"^e_j , 

€2 = sin (p Sji- cos (pSy = |^e""*e+-e"*e_j, (76) 



7 MEAN VALUES OF THE YORP TORQUE PROJECTIONS 

The YORP torques belong to the class of very weak, dissipative factors that may produce significant effects only through their mean values 
that act systematically. If we exclude any commensurability between and a>, we can identify the averaging with respect to time with 

^^^■>=irr^^'^'^' '''' 

rejecting all purely periodic terms of (p and t// in 

Mj = M-ej, 1,2,3. (78) 



7.1 Spin component M3 

According to Eqs. l |72t , the mean component (M3) = M plays the principal role in the YORP effect, being responsible for the systematic 
change in rotation period. Besides, it also can be described with a simpler set of formulae than Mj or M2. For these reasons we consider M3 
prior to the remaining components. 



7.1.1 First order 



The first order terms in M3, to be labeled m9\ result from Af*^'' defined in Eq. ( I66t . But 



(M 



(1)\ 



' ^ 

i>l m=-l 



\ l,ml ' 



(79) 



because accor ding to Eq. (| 54l. (fj'^) = for all m + 0, whereas for m = 0, we have mlj^^^ = 0. Thus, as demonstrated by Nesvorny and 
Vokrouhlicky (|2to3;|2008^ the spin-up or spin-down of an asteroid is a second order effect, involving the products of the shape coefficients. 



7.1.2 Second order 

The second order part of M3 will consist of terms derived from M*''^ and M*"^^ defined in Eqs. l |69l l and i70t . where we substitute from 

Eq. l l54t and f^^ from Eq. l l55t . We begin with M^"' and after the substitution of temperature coefficients J55t with k = m, and g = 0, we 
readily obtain the mean value 

■ 3 ^ w^'^ Hc\ 

(M("))=(M(").,3)=-|!^2 2 /,:.-(4l!)=- z (-i)"-/c«/.>^^o,o4,,^,:r''<o^-'^°"- (80) 

l>l m=-l l,p,j,m 

where 



' corrected 
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3 Vc 



(81) 



1,0 



Recalling that (/jq.o = 0, Kq q = 1, Hq = 2, and w]'^^ = when p is odd, we reduce ^M^"'^ to 



l,J,p^ 1 m=-/ 

using the special case of the overlap i ntegral W„ defined in Eg. ^A27|l. 

In agreement with the results of iNesvornv & Vokrouhlick\i we find that the remaining term M^"^' has no systematic effect on 

the spin rate. Indeed, 

; 



l,j,p,m,q lj,p^l m=-l 



(83) 



thanks to the symmetries between the terms with the indices (l,j,m) = {qi,q2,±q^) and {q2,qi,±q^). For similar reasons we can actually 
replace g[^^ . in Eq. ^ by gf, 

Thus, the complete mean YORP torque (Mj) is given by 



l,j,p^lm=-l 



2p 
0" 



(84) 



Equation iS2l can be converted into an explicitly real form if we combine the positive and negative values of the summation index m 
and use the parity properties of shape coefficients, Wigner d-functions, and Gaunt coefficients. The result, with vanishing terms dropped and 
with optimized summation with respect to I and j, becomes 

; i+j 
{Mi) = -2ff J] J] J](-l)'«mJ(2/ + 2j+ 1) (c,,,„5,+2;,«,-5,,,„Q+2;,m)X V^P^y"'^^^^^^^ 

1 /> 1 m- 1 p= j 

Two remarkable proper ties known from previous works i Capek & Vokrouhlickvll2004 ; Nesvomv & Vokrouhlickvl[2007 ; Scheeresl2007 ; 
Scheeres & Mirrahimil2008 ) are worth noting: (M3) is an even function of (;r/2-e), and it is independent on thermal conductivity. The former 



(85) 



property is a direct consequence of the f act, that <j^^ is actually the Legendr e polynomial P2p(cose). A good argument for the independence 
of (M3) on conductivity can be found in Nesvomv & Vokrouhlickv 1 2008bl) . 



7.2 Obliquity component Mi 

7.2.1 First order 

Using Eq. l l66t . we project M^'-* on e\ defined by Eq. l l76t in terms of e±. Then, recalling that e± - So = e+ ■€+ = S- -S- = 0, and e+ -e- = 1, 
we find that the average value of the first order terms is given by 

l>l m=-l 

Substitution of the spherical temperature function from Eq. J54b . followed by the rejection of purely periodic terms (in both components it 
amounts to setting m = k = 0), leads to 

= -/^_;^/,;;,w;f/fo(/fi,o<T7o4,oe-'>'° + ^^_i,o<orf-i,o^"''^-'°)- (87) 

Using the particular values of the subscripted symbols, namely 



V-l,0 = -¥'l, *^-l,0 = ^l, ^0 = 2, d.'-ifi = -d'xfi^ '^fto = V ^ 9 ' fw-^lfi' 
we find the first order contribution in an explicitly real form 



- — Ki sin^i V Co y/l(I+l)wYU[ „, (88) 

depending only on zonal harmonics coefficients C/q. Moreover, the properties of overlap integrals l lA24t imply that only the even degree 
zonals remain in the final formula, so - introducing the special function Wi from Eq. iA27l - we finally obtain 

(m<'>) = -2aKi sin^pi ^Czz.o V2/(2/+ 1) Wz^O- (89) 
l>l 
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This expressi on confirms the pre sence of the YORP effect in obliquity even for such regular objects as spheroids, provided they have nonzero 
conductivity dBreiter et al ]|2007h . 



7.2.2 Second order 

In the second order we first consider M^^^\ According to Eq. ( I69t . the departure point is similar to Mj'\ i.e., 

m ^ Z Z K ^-") - '^")) • (90) 

l-»\m=-l 

although the temperature coefficient from Eq. l |55t is more complicated. The and ij/ independent part, resulting from k = m and q = 0, 
is 

(M(")) = ,V2a^:, Y,i-irfiJj,rApj[<.ne^^^^^^^^^ (91) 

l,i,pjn 

The treatment of M*"^^ is comparable to MJ ' . The result becomes quite similar to Mj"', because in this case Eqs. l l70t and l l54t lead to 



3v, 

l,J,p,m,q 



Zj J'''»JP''l['^l,ni^l,p.j y j-m-q-lj^'^UnMLpJ \^ j-m-q+ljj- 



il^aK, ^-'rfUrnf;,„\^U&Tp^j'"''e-'^^ (92) 

Z,;,/),m 



In order to add Eq. ([92} to l |91t , we manipulate the former, interchanging the indices {Up,]) , permuting the indices in the Gaunt coefficients, 
and making use of = (-1)"'/^*^^, and cr^_,„ = cr^,„- The resulting total second order torque is 



[Mf^) = (M' 



l,j,p,m 



The marker yu = 1 introduced in Eq. | |7U helps us to trace the influence of the M^^^\ 

After a tedious procedure of combining the terms with (l,j,m) = (qi,q2,±q^) and (l,j,m) = {q2,q\,±q-i), we obtain the explicitly real 

form 

(Mf>) = «^:icos^,;^2^(-l)'"(Q,,„5,,2,„,-Q,2,,„5,„,^^ 

/>1 y>l m=l /j= y 

^ r VVVr iv«Cr 4- v ? \ ^' 2g^p,/+2/2p,/ ~ "^P^^p + l)/z ^^^q ,2p,„., 
+0-^1 sm^pi2^2jZj("^^ [Ci,mCuum+Si„,Su2j,„)—— ^ ^;,/+2/2/, ^1,(^-94) 



Using the recurrence relation l lA32t we were able to simplify the sum of two Gaunt integrals, but it did not help for their difference, which 
explains asymmetry of the two parts of Eq. J94b . Notice, that yu is present only in the terms factored by sim/Ji , i.e., in the part that does vanish 
in the Rubincam's approximation. 

7.3 Precession component M2 

Probably, the YORP induced precession will meet rather limited interest, being negligible when compared to the precession resulting from 
much stronger gravitational torques. Fortunately, the similarity between Si and 62 resulted in a very simple rule: it suffices to apply a 
formal substitution (s'mip[, cos(/3i) —> (-cosipi, s'mip[) to turn m|'' or m|^^ into M2''' or M^^ respectively. Of course, this observation is 



based on the complete derivation, that we omit for the sake of brevity, but it can also be deduced from the equations of lScheeresI ( 120071) or 



Scheeres & Mirrahimi. C2008.) . 



7.4 Summary 

In order to collect the final equations for the average YORP torque in a possibly compact form, revealing their actual structure, we first 
introduce few auxiliary symbols 

^ i.m ~ ^l,mS 1+2 j,m ~ S l^mCl+2 j,m^ (95) 
^l.m ~ CijnCl+2j,m+ S ijnS i+2j,m^ 

(96) 

and 
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si=Kis'mipi, ci=Kicosipi, 



(97) 



that should not be confused with c = cose, and s = sins. 

Further, we note that only the special cases of Wigner d-functions are required: rfy y and rf^'^j, which are simply related with the Legendre 
polynomials P2p(c) and Legendre functions f ^^^(c), respectively. We will also limit the summation according to a practical assumption, that 
shape harmonics are known only up to some finite degree and order A'. Thus, using Eq. jA22t . we conclude this section with the following 
set of formulae: 



(Ml) 
(M2) 



-la 



l+j 



2 XI 2 ^im 2 



/=! j={) m=0 



p=max(l,;) 



l=\ j=\ m=l p=j 



^-24^] / i+i 

"Y. Z Z^i«Z^'."W.P^2;Xc), 
/=! >1 m=l 



{Mi) 



wher^ 

l.m,j,p 

l,m,j,p 
^l,m,j,p 



P=J 



-2(-ir jm(2/+2;+ w, 



(-1)' 



2-5o,„ 



(J 



^IpJ+ljSipJ 



l,l+2j,2p "P' 
l+2j 



p{2p+\) 



-2p 



^lM2j,2p^P' 



(-l)'«7(2/ + 2;+l) 



V2p(2p+1) 



'Ml,m,j,p M'/j, 



(98) 
(99) 

(100) 

(101) 

(102) 
(103) 



is the entier function, rounding down to the nearest integer, and Sj^i^ is the Kronecker delta function. 

The coefficients Wp are quickly decreasing in absolute value from Wq = - ■0r/4, with limp^oo Wp = 0. They can be evaluated by means 
of the simple recurrence 



W„ 



4p+l (2p- l)(2p-3) 



Wp-i. 



[4p-3 2p(2p + 2) 

The Gaunt coefficients present in the final expressions are defined in terms of the Wigner 3-j symbols in Eq. l |A29t , so 



>->m,-m,0 
^l,l+2j,2p ' 



and 



(2/+l)(2/ + 4;+l)(4p+l) 



4n 



I l + 2j 2p 




I l + 2j 2p 
m -m 



l,mj,p ■ 



(2Z+l)(2/ + 4;+l)(4p+l) 



2;r 



l + 2j 




2p 




/,m 



/ 

m— 1 



l + 2j 



2p 
1 



I 

m + I 



l + 2j 



2p 
-1 



(104) 



(105) 



(106) 



The 3-j symbols and Gaun t coefficien ts admit a n umbe r of recurrence relation and their computation is a routine task in physics (eg. 
Luscombe & LubanI ( Il998h . Ixiil Jl996l) , Isebilleaul ( Il998l) ). A reliable software for the 3-j coefficients can be found in the open source 
SHTOOL^ package, although major computer algebra systems implement the 3-j symbols as one of standard functions. The values of 



■m,~m,0 
l+2j,2p 



and "J^imjp decrease with the growing p, but they do not reveal any monotonicity with respect to /, n, or j. 



The presence of even degree Legendre polynomials of c in (M3) and of Legendre functions P\p{c) in (Mi_2)(fi) explains the symme- 



tries well known from dRubincamlbood) or JCapek & Vokrouhlickil2004l) : (MiXs) = {Mi){n-E) and {Mi){e) = -{Mi){n-s). Similarly, 
{M2){e) = -{M2)(n - e), although the division of M2 by s in the right hand sides of Eqs. l l74t and J75t reverts this symmetry and the effect 
on precession is an even function of the angle (| - e). 

In the Rubincam's approximation, when we negl ect the conductivity, i-i = and ci =1. This simplification h as no effect on the 
spin component {Mo,}, in agreeme nt with the results of ICapek & Vokrouhlickv (2004); Scheeres & Mirrahimi 1 2008h : Mvsen 1 2008h and 
Nesvomv & Vokrouhlick'^ J2008bl) . 



* corrected 

^ available at http://www.ipgp.jussieu.fr/~wieczor/SHTOOLS/SHTOOLS.html 
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8 LEGENDRE SERIES FORM AND TRUNCATION RULES 



If the Legendre functions in the formulae of Sectionl7.4l are expressed in terms of the powers of c, we obtain exactl y the same final expressions 
for sp in component (M3) as lNesvornv & Vokrouhl ickv (2007). Inspecting the coefficients of such expansions in jNesvornv & Vokrouhlickv 
20071) one finds no rule that might help to select the leading terms, responsible for the major part of the YORP torque. The situation is much 
better in the present formulation, because the monotonicity properties of Wp alone (even fortified in the products S'J'l^^j'lp '^lj,i".p ^p) 
suggests to attach the major importance to the terms with the lowest p values. In this c ontext, we can justify the partition of the YORP 
torque terms into subsequent 'orders' of j, proposed bv lNesvomv & Vokrouhlick^ 1 2007[), because th e index j establishes the lower bound 
on the range of p. Truncating with respect to j is also consistent with the approach of MvsenI (2008), who imposed it at the beginning of 
his derivation, using an insolation function that actually is restricted to the first two Legendre polynomials Pi and P2, with the efi'ect of Pi 
disappearing after the averaging. 

The sums presented in the previous section can be easily rearranged from the form attaching to each shape coefficients combination a 
sum of Legendre functions, to an alternative form of Legendre series with a sum of coefficients attached to each Legendre function. Simple 
manipulations lead to 



(Ml) 



{M2} 



{Mi) 



Y^[siX\ + SlXl+ClZq) P\^{C), 
N 

Y^{-ciX\-ciXl + siZq) P\^(c), 



q=l 
N-1 



0'YjAgP2q{c), 

q=l 



with the coefficients given as the sums 



q N-2j I ^ q N-lj I 

^? " Zi Zi Zi '^im ^ImJ.q^ = Z Z Z 



j=0 1=1, m={) 

where the lower summation bound is Ii 



,/=I 1=1, m=\ 



q N-2j I 

= Z Z Z ^im ^''"'•J'l 
j=l 1=1, m=l 



(107) 



(108) 



(109) 



(110) 



■ max(l, q-j). 



These Legendre series suggest to introduce a notion of the 'YORP degree' , complementary to the 'YORP order' proposed bv lNesvomv & Vokrouhlickv 



( 20070 . A term of the YORP degree ^ is a sum of terms with YORP orders j < q. Truncating the YORP series one should first decide on the 
maximum degree ^max, and then either maintain the maximum order j naturally bounded by the current value of q, or additionally impose 
some smaller value of jmax- According to our experience, a decent approximation level is reached already at qm^x = 2, so a qualitative pattern 
of the E dependence is established in the second degree truncated series^ 



9 TEST APPLICATIONS 



9.1 General remarks 

The basic assumptions of our solution are similar to the ones of iNesvomv & Vokrouhlicki ( l2007Ll2008bl) : almost spherical homogeneous 
bodies with constant albedo and small thermal penetration depth combined with low conductivity. The shape restri ction is quite severe, ruling 
out m any irregularly shaped objects like (433) Er os. A good account of this difiicul ty was given in Section 8 of ( Nesvornv & Vokrouhlickvl 
2007). Our experience confirms the statement of iNesvomv & Vokrouhlickv J2007h that qualitatively correct results can be obtained even 
for such irregular bodies, although the main cause of discrepancies is not the magnitude of the shape coefficients as such, but rather the 
occurrence of major concavities in the figure of a considered body that generates large shadow zones treated as illuminated in analytical 
insolation function models. So, we test our solution on a spheroid, that is obviously convex, and on the asteroid 1998 KY26 - a reasonably 
regular object. 



9.2 Spheroids 



Breiter et al 



Breiter et al 



As a first test we confronted the present results with YORP torque formulae for spheroids provided by 
results for spheroids used a simplistic thermal lag model of 'delayed Sun', the final formulae from Sect. 3.5 of 
easily generalized to a more realistic thermal model by replacing sind with i i of the present paper Then, with the shape model 

C4,0 ' 



( 20070 . Altho ugh the 
( I2OO7I) can be 



2 11 / 
e + — e 
21 



2-v9r 4 



(111) 



* Coefficients files up to W = 50 are available from the first author (SB) upon the request. 
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Figure 2. Asteroid 1998 KY26 seen from the north pole with a simulated lighting from the equator. The regions where our insolation function model fails are 
labeled A (self-shadowing) and B (terminator cuts olf a sunlit region). 



derived for an oblate spheroid 



ae Vl-e2 
Vl-e2w2' 



with the semi-axis a^. related to the mean radius a by the series in eccentricity e = -yJT^-cla^ 

a a; flp i 

\ 6 120 

we find from Eq. (ToTJ 



On the other hand, according to Breiter et alj j2007 ) 



71 (l-A)(S>al loll 5,14 



(112) 



(113) 



(114) 



(115) 



Replacing s'm6 by si and substituting Eq. ill3t we obtain a perfect agreement between Eqs. i ll 141 ) and l lll5t . 

The comparison is meaningful only up to e'^, because = Oie"^) marks the limit of our present second order approximation. For the 
same reasons we dropped C^^o = 0(e^) and higher degree zonal coefficients. 



9.3 Asteroid 1998 KY26 



From the collection of numeri cally simulated mea n YORP torques presented b y ICapek & Vokrouhlick\1 ( |2004|) . we selected the asteroid 
1998 KY?fi, c onsidered also by Scheeres 1 2007 ) and Scheeres & Mirrahimi 1 200 8h . The same object was used by Nesvorny and Vokrouhlicky 
2007|j2008bl as an example for their second order model. The authors were quite optimistic about this body, and indeed, a side view presented 
in their paper looks like well fitting the assumption of the theory. However, the polar view (Fig. ^ reveals a source of possible problems 
when the obliquity e is close to or n. Some regions, qualified as sunlit by the spherical terminator criterion are actually in shadow, whereas 
some sunlit regions are neglected and considered dark. Apparently, these two effects cancel out in the numerical simulation, but the analytical 
model tends to sharpen their influence. 

Nevertheless, we constructed the spherical ha rmonics model up to degree and order N = 26 using the least squares adjustment on the grid 
of 2048 surface points determined bv lOstro et al.l ( 11999;) . reduced to the reference frame of centre of mass and principal axes of inertia. The 
standard deviation of our harmonics set on the grid points was 0.2 mm, quite satisfactory for an object with the mean radius a = 13.1913 m. 
Interested rather in c omparison than in the dy namics of the asteroid itself, we used the same, partially fictitious values of physical parameters 
that were applied by Capek & Vokrouhlicky ( 2004) : rotation period 6 hr, density p = 2.5 g/cnr', orbital semi-axis Oq = 2.5 au, heat capacity 
Cp = 680 J/(kgK), and albedo A = 0. The maximum moment of inertia, that we derived during the reduction to the principal axes, was 
J3 = 1.7273 X lO^kgm^. 
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Figure 3. Present solution (continuous lines) for 1998 KY26 and the numerical model o f lCapek & Vokrouhlickvl J2004 (dashed lines). Thick lines refer to the 
conductivity K = 0. Thin lines present the results for K = 0.001 and K = 0.01 W/(mK). 



The agreement of 6j computed from Eq. ( |109b shown in Fig. [3] with the resuUs of ^Capek & Vokro uhlickv 1 2004*) is as good as the one 
obtained by Nesvorn y & Vokrouhlicky. (2007). It is no t surprising, because we expanded our expressions in power series of c and compared 
them with Table I of iNesyornv & Vokrouhlicky! <2007l) . concluding a perfect agreement of the nu merical coefficients. Yet, when we pass to 
the obliquity yariations e, we find much more distortion with respect to smooth curyes taken from lCapek & Vokrouhlickvl (|2004|) . although 
the orders of magnitude and qualitatiye dependence on obliquity and thermal conductivity may be qualified as correct. In our opinion, 
the degraded accuracy at small values of £ or (tt - e) should be attributed to the north-south directed ele v ated p atterns in the topography 
of 1998 KY26- We have observed, that the departure from the sine-like curve of ICapek & Vokrouhlickvl | |2004|) strongly depends on the 
maximum degree and order N of the accounted shape harmonics. As a matter fact, a sufficiently inaccurate shape model, with the elevated 
regions smoothed out, gives better shapes of the s curve, albeit producing a smaller amplitude. 



10 CONCLUSIONS 

Although we used the same assumptions and general strategy as Nesvorny and Vokrouhlicky (2007; 2008bl) , we have obtained a signifi- 
cantly simpler form of the final expressions of the mean YORP torque. In particular, the final formulae are expr essed in terms of Legendre 
polyno mials and associated functions, as compared to the products of Wigner functions that can be found in iNesvomv & Vokrouhlicky 
1 2007 ). The advantage comes from the way we expand the insolation function ~ more straightforward than the one exposed in Appendix B of 



Nesvorny & Vokrouhlickvl (120070 . Our formulation is well suited for the partition into the terms of primary and secondary influence accord- 



ing to the YORP degree, directly related with the degree of Legendre functions present in the final expressions. Moreover, it has significantly 
better numerical properties than the power serie s of c. We confirm the basic facts established in recent years, lik e the independence of the 
spin related torque component on conductivity I' Capek & Vokrouhlicky 2004 ; Nesvorny & Vokrouhlick y '2008b') within the plane parallel 
temperature model, or the presence of YORP elfect in obliquity for the bodies with rotational symmetry ( Breiter et al. 2007). The explicit 
appearance of Piic) in Eq. ( I109t and the remarks about the properties of Wp, neatly explain the fact that the spin related YORP torques tend 
to vanish close to the zeros of /'2(c) dNesvornv & Vokrouhlickvll2007D . Similarly we can predict, that if the obliquity component (Mi) has 
zeros different than c = or ^ = (implied by each Legendre function P^^^ic)), the new zeros should be located in the domain 40° < e < 60° 
(or 120° < e < 140°) implied by P\{c) ± P\{c) = 0. This conjecture is confirmed by a large sam ple of G aussi an random spheres s imulated in 
I Capek & Vokrouhlicky 12004 ). as well as by the observed 'Slivan states' in the Koronis family ( Slivanlliooi : Vokrouhlicky et al ]|2003h . 

The assumptions about the insolation model seem to be the main sourc e of divergence between numerical models and a completely 
analytical formulation. In this context, the approach taken by iMysen J2OO8I) seems to be a good compromise: in his paper the insolation 
function is derived from numerical simulation and then incorporated into an analytical treatment. As usually, however, such semi-analytical 
approach breaks the direct link between the shape coefficients and the YORP torque magnitude. It is a hard duty of an analytical solution to 
show such link, even if the final relation becomes as cumbersome as the one derived here. 

Notably, the YORP effect in the Rubincam's approximation is a second order elfect in terms of the shape coefficients. Which of the 
harmonics products will play the leading role, depends on the particular body shape, because there is no monotonicity of the numerical 
coefficients accompanying a particular harmonics product. It is easier to point out nonsignificant harmonics: for example, the spin torque 
(M3) is independent on the zonal harmonics coefficients C/.o. Although the presence of the nonzero conductivity brings in the first-order 
terms to the (Mi) and (M2), in most cases they will be less significant than the second order contribution, because si - factoring the zonal 
coefficients - is typically much smaller than C]. Thus the zonal harmonics coefficients are generally of minor importance, unless the body 
has a highly regular shape, close to a rotationally invariant solid. 
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The modular structure of the present work allows a modification of particular elements of the solution without restarting the entire 
work from the scratch. Our next goal is the application of a more advanced temperature solution. The results are promising and will soon be 
announced in a separate paper. 
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APPENDIX A: SPECIAL FUNCTIONS 



Generally speaking, we use the definitions of special functions according tol Biedenham & Louckl ( 119811) . The only exception from this rule 
is the inclusion of the Condon-Shortley phase (-1)"' in the definition of associated Legendre functions P'". All formulae from other sources 
have been adjusted to fit this convention. 

Al Spherical harmonics 

If Pi(u) stands for the usual Legendre polynamial of degree /, we define the associated Legendre functions (ALE) as 
Introducing the normalization factor tr/ ,„ 



21+1 {l-m)\ 

o-l,m = A — ; — 77" — -7, (A2) 
V 4;r (l + my. 
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we define the normalized ALF 
Q",'{u) = a,^„Pl\u). 

Spherical harmonics 7/ „, used in this paper are 

where u is the sine of latitude (i.e., cosine of colatitude) and A is the longitude of the point on a unit sphere. 
The symmetries of Legendre functions 

(l-mV 

' (t + my. ' 

@j'"(u) = (-l)'"@1'(u), 

lead to 

y;,„(M,^) = (-irr,,-„,(M,i), 

where an asterisk indicates a complex conjugate. 
The following orthogonality relations hold 

X2n n 1 r^ln r\ 

AX J Y j^i,{u,X)YlJu,X)Au = (-\T AX J Yj^kiu.X)Yi^_„(u,X)Au = Sji8i,^,n, 

with Sp^q designating the Kronecker delta function equal \o\'\i p = q and zero otherwise. 

A function X) defined on a unit sphere can be expanded in spherical harmonic series 



(A3) 
(A4) 

(A5) 
(A6) 

(A7) 
(A8) 



F{u,X) = 22 a;,,„K/,„,(«,i), 
/>0 m=-/ 

dM J F{u,X)Y*iJyU,X)diX. 
If F is real-valued, then the coefficients a/ „, admit the property similar to l |A7t 

"/m = (-l)"'"/,-m- 



(A9) 
(AlO) 

(All) 



A2 Wigner functions (generalized splierical liarmonics) 

Consider two Cartesian reference frames Oxyz and Ox'y'z' having the same origin O, and related by rotation defined in terms of the 3-2-3 
Euler angles a, /3, y 



r' = R3(7)R20e)R3(ff)r, 
with the standard rotation matrices 
(10 



Rl(0): 



cos(^ 
-sini! 



sin^ 
cost 



R2(^) = 



cos^ 




-sin^ 

1 



R3(<^) ■■ 



cos cj) sincj} 
-sin0 COS0 
1 J 



sin^ cos^ 

The same transformation can also be specified in terms of more common 3-1-3 Euler angles 
r' = R3(7-;r/2)Ri(/?)R3(a + ;r/2)»-. 



(A 12) 



(A13) 



(A14) 



According to the transformation laws of spherical harmonics for the 3-2-3 Euler angles sequence, each spherical function of polar 
coordinates (d,X) in Oxyz frame becomes a combination of the same degree spherical functions of polar coordinates 9' and X' in Ox'y'z'- 
Thus, for each harmonic 



y,,„(cose', X') = 2 Dlja,p,j)Yi^kicose, XI 

k=-l 

or, conversely, 

y,,;„(cos0, X)=Y, [^L,t(^.Ar)]' n-tCcose', x'). 



(A15) 



(A 16) 



k=-l 



^ corrected 



1 8 Breiter and Michalska 



where the general definition of Wigner D-matrix elements (or Wigner D-functions) is i Biedenham & Louc3ll981 ) 

d',„/<^A7) = dl,(P)e-'^'""^'^\ (A17) 
The symbols j^ifi) designate so called Wigner d-functions 



<*03) = A TT^^^tItJ-^ [smOS/2)]*-'"[cosC8/2)]*-« pf-^"'''^"'\cos/i). (A18) 
V (/ + m)!(/-m)! ' * 

The definition of d-functions in Eq. JAlSt involves the Jacobi polynomials 

'rv.=z(:r)(:;::)(iiirm"' 



Wigner d-functions ^ admit some parity properties like 



^'. = (-1)'""^' rfL = (-l^'"4..- (A20) 



Using these two identities, one can always reduce the ^ function to the case when k + m^Q and k-m^Q, avoiding the appearance of 
negative powers of sines and cosines in Eq. JAlSb . 

As it follows from ( IA20t . the conjugate of ^ is 

[Dl ,{a,p,y)] = (-1)'«-*Z)'_„, 7). (A21) 

In the special case rf^ ^ or dy ^, Wigner d-functions reduce to AFL 



A3 Overlap integrals and Gaunt coefficients 

The Wigner 3-j symbols are numerical coefficients defined as a sum 

(jl +J2-73)!0'l - J2+ j3)!(-;i +,/2 + ,/3)! 



Jl J2 J3 
mi m2 mj, 



(jl+72 + ;3 + l)! 



VO'l + "II ) ! ( jl - mi ) ! 0'2 + ^2) ! O2 - "J2)! (73 +m-i)\ {j3 -1113)1 



k^- (71 + 72 - 73 - ky. (71 -mi-ky. (72 + OT2 - fe)! (73 - 72 + m + ky. (73 - 71 " '"2 + fc)! 



(A23) 



where the index k takes all integer values leading to nonnegative factorial arguments. 

The so-called overlap integral is an integral of the product of two associated Legendre functions with the integration limits ±1. We use 
a normalized version of this function 

<';r=£c;(")0;?(«)dM. (A24) 

Using a formula of Dong & LemusI 1 2002h . converted to the normalized version, we can express the overlap integrals in terms of the Wigner 
3-j symbols 



(-1)' 



min(m 1 ,m2) 



2l'«i-'"2|-4 |„j _ V(2ni + l)(2;J2 + l) 



^ (2/+l)G(|mi-m2l,0 
/=max(|«i -«2l.k"i ^'"2!) 



ni 112 I 




ni 112 I 
—m\ 1112 m\-m2 



G{m,D = (l+(-l)'+"') 



(/-m)! r(//2)r((Z + m+l)/2) 
(/ + m)! ((/-m)/2)!r((/-l-3)/2)' 



(A25) 



(A26) 



Obviously, the overlap integrals are symmetric with respect to the indices w^'jj^^ = w^^jj^', and they vanish either if |m/| > n,-, or if the 
integrand in Eq. ( IA24l l is an odd function of u, i.e., when (-l)"i+«2+"ii+ni2 - _i Moreover, because of the orthogonality of the associated 
Legendre functions, we have w™'™, = 0, for nj 1 112- 

Because of a frequent appearance of w'j'j^^ cr in the final formulae we introduce the special case of the overlap integral 



'^11 V^ / (2p-l)!! 
0-1,1 4''+! \ p\ 



V^(4p+1) 

(p+l)(2p-l) " 4(p-H)(2p-l) 



\P2pm\ 



(All) 



The Gaunt coefficients (or Gaunt integrals) are defined as an integral of three spherical functions over the unit sphere surface 



Their expression in terms of the Wigner 3-j symbols is much simpler than for the overlap integrals 
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(A28) 



■2m\,m-,,mi, 

1,(2,(3 



(2/1 + 1x2/2 + 1x2/3 + 1) 



An 



h h h 





h h h 
mi m2 m-i 

Gaunt coefficients are non-zero if the following three conditions are satisfied: 
(-1)''+'=+'' = 1, mi+m2+m3=0, \li - lj\ ^ k ^ li + Ij. 

Note, that the signs of all m,- superscripts can be simultaneously inverted, and the pairs /, ,«?,■ can be arbitrarily permuted 



l2,h ~"h,h,h ~"l2.h,l} '^^hJuk 

From the three-term recurrence identities related with the Wigner 3-j symbols i Luscombe & LubMilll998 ), we can deduce useful relation 
(/3(/3 + 1) - /2 (/2 + 1) - /l (h + D) S'fXi'^ = y/h(h + i) [ ^|{h-m)(h+m+\)gl'-^l^-^^'-'" + + m^i - m+ Dff'^Jf^^'-"'] . (A32) 

The Gaunt coefficients appear in the transformation rule for the product of spherical harmonics (Clebsch-Gordan series). Applying Eq. ( IA28t 
to l lAlOt we find 



^ etc. 



(A29) 



(A30) 



(A31) 



Ji+J2 



Yj,,kjj2.k2= 2j 

/=max(|>i-,/2l,l*i+/:2l) 

where the sum ji + ji + l must be even. 



(-1) ,,/,,/ 'lM+k2 



(A33) 



